//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <dirent.h>
#include <math.h>

#include <iostream>
#include <fstream>

#include "zlib.h"

#include "libraries.h"
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
#define LINGOLENGTH  4
#define LENGTHARRAY  40096
#define CUTOFFFREQ   1
#define ORIGINALFREQ 0
#define CHANGEDIGIT  '?'
#define EUCLIDIAND   1
#define POSITIONL    1
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
using namespace std;
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
struct SmilesRef {

  String* smiles;
  String* ref;

  SmilesRef() {
    this->smiles = new String();
    this->ref    = new String();
  }

  SmilesRef(char* smiles, char* ref) {
    if (smiles && ref) {
      this->smiles = new String(smiles);
      this->ref    = new String(ref);
      delete[] smiles;
      delete[] ref;
    }
  }

  ~SmilesRef() {
    if (this->smiles) delete this->smiles;
    if (this->ref)    delete this->ref;
  }

  struct SmilesRef operator=(const struct SmilesRef &smilesRef) {
    this->smiles->In(smilesRef.smiles);
    this->ref->In(smilesRef.smiles);
    return *this;
  };

};

typedef struct SmilesRef SMILESRef;
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
double obtainLingoSim(char *lingoString1, char *lingoString2) {

  if (strlen(lingoString1) < LINGOLENGTH || strlen(lingoString2) < LINGOLENGTH) return 0;

  int length1 = strlen(lingoString1), length2 = strlen(lingoString2);
  char tmpLingo1[length1], tmpLingo2[length2];
  int  i, j, k, z;

  //Change the Cl bu L and Br by R in the two lingo sequences
  for (i = 0, j = 0, k = 0; i < ((length1 > length2)?length1:length2); i++) {
    if (i < length1) {
      if      ((lingoString1[i] == 'l') || (lingoString1[i] == 'r')) tmpLingo1[j - 1] = (lingoString1[i] == 'l')?'L':'R';
      else if (isdigit(lingoString1[i])) tmpLingo1[j++] = '?';
      else tmpLingo1[j++] = lingoString1[i];
    }
    if (i < length2) {
      if ((lingoString2[i] == 'l') || (lingoString2[i] == 'r')) tmpLingo2[k - 1] = (lingoString2[i] == 'l')?'L':'R';
      else if (isdigit(lingoString2[i])) tmpLingo2[k++] = '?';
      else tmpLingo2[k++] = lingoString2[i];
    }
  }
  tmpLingo1[j] = '\0'; tmpLingo2[k] = '\0';

  j -= (LINGOLENGTH - 1); k -= (LINGOLENGTH - 1);
  int    L = j + k; //Total Lingos in two sequences
  double similarityValue = 0, Nai, Nbi;

  for (i = 0; i < j + k; i++) {

    char   compare[5];
    bool   analize = true;

    Nai = 0, Nbi = 0;

    if (i < j) {

      strncpy(compare, tmpLingo1 + i, 4);
      compare[4] = '\0';

      for (z = 0; z < ((j > k)?j:k); z++) {
        if (z < j && !strncmp(compare, tmpLingo1 + z, 4)) {
          if (z < i) {
            analize = false;
            break;
          }
          Nai++;
        }
        if (z < k && !strncmp(compare, tmpLingo2 + z, 4)) {
          Nbi++;
        }
      }
    }
    else {

      strncpy(compare, tmpLingo2 + i - j, 4);
      compare[4] = '\0';

      for (z = 0; z < ((j > k)?j:k); z++) {
        if (z < j && !strncmp(compare, tmpLingo1 + z, 4)) {
          analize = false;
          break;
        }
        if (z < k && !strncmp(compare, tmpLingo2 + z, 4)) {
          if (z < i - j) {
            analize = false;
            break;
          }
          Nbi++;
        }
      }
    }

    if (analize == true) similarityValue += (((Nai - Nbi) >= 0)?(1 - ((Nai - Nbi)/(Nai + Nbi))):(1 - ((Nbi - Nai)/(Nai + Nbi))));
    else                 L--;
  }

  return similarityValue/L;

}
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
void trim(char *stream) {

  int i, j;
  for (i = 0, j = 0; i < strlen(stream); i++)
    if ((stream[i] != '\n') && (stream[i] != '\r'))
        stream[j++] = stream[i];
  stream[j] = '\0';

}
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
char* splitPosLine(char *line, int pos, char character) {

  if (pos < 0) return NULL;

  int   count = 1, length = 0;
  char *tmpResult = new char[strlen(line) + 1];

  for (int i = 0; i < strlen(line); i++) {
    if (line[i] == character) { count++; i++; }
    if (pos == count) tmpResult[length++] = line[i];
  }
  tmpResult[length] = '\0';

  char *result = new char[length + 1];
  strcpy(result, tmpResult);
  result[length] = '\0';
  delete[] tmpResult;

  return result;

}
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
int main(int argc, char *argv[]) {

  String* nameSeqFile = NULL;
  String* nameOutFile = NULL;
  String* location    = NULL;
  int     radius      = 0;

  if (argc < 5) {
    cout << "Error in the input arguments. The execution must be: " << endl;
    cout << "./ComputeAverageNeighborSequences <Input Ordered SMI file> <Input Location Similarity Files> <Output file - GZ> <Average Radius>" << endl;
    return 0;
  }
  else {
    nameSeqFile = new String(argv[1]);
    location    = new String(argv[2]);
    nameOutFile = new String(argv[3]);
    radius      = atoi(argv[4]);
  }

  if (!nameSeqFile || !nameOutFile) {
    if (nameSeqFile)   delete nameSeqFile;
    if (nameOutFile)   delete nameOutFile;
    if (location)      delete location;
    return 0;
  }

  char                buffer[LENGTHARRAY + 1];
  TListE <SMILESRef>* listSequencesRef = new TListE <SMILESRef> ();

  gzFile sequenceFile = gzopen(nameSeqFile->Out(),   "rb9");
  while (!gzeof(sequenceFile)) {
    gzgets(sequenceFile, buffer, LENGTHARRAY); trim(buffer);
    if (!buffer) continue; 
    SMILESRef *tmpSmilesRef = new SMILESRef(splitPosLine(buffer, 1, ' '), splitPosLine(buffer, 2, ' '));
    tmpSmilesRef->smiles->Trim();
    tmpSmilesRef->ref->Trim();
    listSequencesRef->Add(tmpSmilesRef);
  }
  gzclose(sequenceFile);

  double AverageLingoSimAverage = 0;
  gzFile outFile                = gzopen(nameOutFile->Out(),   "wb9");
  
  for (int i = 0; i < listSequencesRef->Length(); i++) {
          
    double lingoSimAverage = 0;
  
    if (i - radius < 0) {
      for (int j = 0; j < i + radius; j++) 
        if (i != j) {
          lingoSimAverage += obtainLingoSim(listSequencesRef->Get(i)->smiles->Out(), listSequencesRef->Get(j)->smiles->Out());
        }
      lingoSimAverage /= i + radius;
    }
    else if (i + radius >= listSequencesRef->Length()) {
      for (int j = i - radius; j < listSequencesRef->Length(); j++) 
        if (i != j) {
          lingoSimAverage += obtainLingoSim(listSequencesRef->Get(i)->smiles->Out(), listSequencesRef->Get(j)->smiles->Out());
        }
      lingoSimAverage /= (listSequencesRef->Length() - i - 1 + radius);
    }
    else {
      for (int j = i - radius; j <= i + radius; j++) 
        if (i != j) {
          lingoSimAverage += obtainLingoSim(listSequencesRef->Get(i)->smiles->Out(), listSequencesRef->Get(j)->smiles->Out());
        }
      lingoSimAverage /= 2*radius;
    }

    gzprintf(outFile, "%s %08f\n", listSequencesRef->Get(i)->ref->Out(), lingoSimAverage);
    printf("%s %08f\n", listSequencesRef->Get(i)->ref->Out(), lingoSimAverage);

    AverageLingoSimAverage += lingoSimAverage;
    
  }

  AverageLingoSimAverage /= listSequencesRef->Length();

  printf("Average LingoSim Average %f\n", AverageLingoSimAverage);
  
  gzclose(outFile);
  
  if (listSequencesRef) delete listSequencesRef;

  if (nameSeqFile)      delete nameSeqFile;
  if (nameOutFile)      delete nameOutFile;

  return 1;

}
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
//---------------------------------------------------------------------------------------------------------------------------------------------------------//
